# Title     : TODO
# Objective : TODO
# Created by: Administrator
# Created on: 2019/6/18

library(ggplot2)
library(dplyr)
library(tidyr)
library(tibble)
library(gridExtra)
library(ggpubr)
library(egg)
library(lemon)
library(optparse)

createWhenNoExist <- function(f) {
  !dir.exists(f) && dir.create(f)
}

option_list <- list(
  make_option("--i", default = "AllMet.csv", type = "character", help = "metabolite data file"),
  make_option("--g", default = "SampleInfo.csv", type = "character", help = "sample group file"),
  make_option("--sc", default = "sample_color.txt", type = "character", help = "sample color file"),
  make_option("--config", default = "config.csv", type = "character", help = "config file")
)
opt <- parse_args(OptionParser(option_list = option_list))

options(digits = 3)

configData <- read.csv(opt$config, header = F) %>%
  set_colnames(c("arg", "value")) %>%
  column_to_rownames("arg")

configData

isPaired <- configData["isPaired", "value"]

sampleInfo <- read.csv(opt$g, header = T, stringsAsFactors = F)
classNotes <- sampleInfo %>%
  .$ClassNote %>%
  unique()
cn <- combn(classNotes, 2)
for (i in 1:ncol(cn)) {
  row <- cn[, i]
  group1Name <- row[1]
  group2Name <- row[2]
  parent <- paste0(group2Name, "_", group1Name)
  createWhenNoExist(parent)
  pcDataFileName <- paste0(parent, "/PCA_", group2Name, "_", group1Name, "_data.txt")
  pc12 <- read.table(pcDataFileName, header = T, sep = "\t", stringsAsFactors = F, comment.char = "") %>%
    as.data.frame() %>%
    select(c("SampleID", "p1", "p2")) %>%
    left_join(sampleInfo, by = c("SampleID"))
  if (isPaired) {
    pc12 <- pc12 %>%
      arrange(Pair_ID)
  }

  fileName <- paste0(parent, "/PCA_", group2Name, "_", group1Name, "_Boxplot_First2PCs.pdf")

  pdf(fileName, width = 10, height = 10)

  parameterFileName <- paste0(parent, "/PCA_", group2Name, "_", group1Name, "_Parameters.csv")
  parameterData <- read.csv(parameterFileName, header = T, stringsAsFactors = F, comment.char = "")

  impoPc1 <- parameterData[1, "PC1"]
  impoPc1 <- round(impoPc1 * 100, 2)
  impoPc2 <- parameterData[1, "PC2"]
  impoPc2 <- round(impoPc2 * 100, 2)

  top <- 0.2
  left <- 0.3
  right <- 0.2

  sampleColDf <- read.table(opt$sc, header = T, sep = "\t", stringsAsFactors = F, comment.char = "") %>%
    select(c("ClassNote", "col"))
  sampleCols <- sampleColDf %>%
    deframe()

  p <- ggplot(pc12, mapping = aes(x = p1, y = p2)) +
    xlab(paste("PC1(", impoPc1, "%)", sep = "")) +
    ylab("") +
    theme_bw(base_size = 8.8, base_family = "Times") +
    theme(axis.text.x = element_text(size = 9, vjust = 0.5), panel.grid.major = element_blank(),
          panel.grid.minor = element_blank(), axis.text.y = element_text(size = 8.8, angle = 90), legend.position = "none",
          axis.title.y = element_text(size = 11), axis.title.x = element_text(size = 12),
          panel.border = element_rect(size = 0.75), plot.margin = unit(c(top, right, 0.5, left), "cm")
    ) +
    #0 line
    geom_vline(aes(xintercept = 0), colour = "#BEBEBE", linetype = "solid") +
    geom_hline(aes(yintercept = 0), colour = "#BEBEBE", linetype = "solid") +
    #point
    stat_ellipse(colour = "#BEBEBE", size = 0.3, level = 0.95, type = "norm") +
    geom_point(aes(colour = factor(ClassNote)), size = 4, stroke = 0) +
    #lengend
    scale_colour_manual("", values = sampleCols)

  rTest <- compare_means(p2 ~ ClassNote, data = pc12, method = "wilcox.test", paired = isPaired)
  rPValue <- rTest$p %>% round(3)

  rBoxPlot <- ggplot(pc12, mapping = aes(x = ClassNote, y = p2, fill = ClassNote)) +
    ylab(paste("PC2(", impoPc2, "%)", sep = "")) +
    theme_bw(base_size = 8.8, base_family = "Times") +
    theme(axis.text.y = element_text(size = 8.8), legend.position = 'none', axis.title.y = element_text(size = 11),
          legend.margin = margin(t = 0.3, b = 0.1, unit = 'cm'), legend.text = element_text(size = 6),
          axis.title.x = element_text(size = 12), panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
          axis.ticks.x = element_blank(), panel.border = element_rect(size = 0.75),
          plot.margin = unit(c(top, right, 0.5, left), "cm")
    ) +
    geom_boxplot() +
    #lengend
    scale_fill_manual("", values = sampleCols) +
    scale_x_discrete(paste0("Pvalue=", rPValue), labels = rep("", length(pc12$ClassNote)), breaks = pc12$ClassNote)

  bTest <- compare_means(p1 ~ ClassNote, data = pc12, method = "wilcox.test", paired = isPaired)
  bPValue <- bTest$p %>% round(3)

  bBoxPlot <- ggplot(pc12, mapping = aes(x = ClassNote, y = p1, fill = ClassNote)) +
    xlab("") +
    ylab("") +
    theme_bw(base_size = 8.8, base_family = "Times") +
    theme(axis.text.x = element_text(size = 9, vjust = 0.5), panel.border = element_rect(size = 0.75),
          axis.text.y = element_blank(), legend.position = 'none', axis.ticks.y = element_blank(),
          axis.title.y = element_text(size = 12), legend.margin = margin(t = 0.3, b = 0.1, unit = 'cm'),
          legend.text = element_text(size = 6), axis.title.x = element_text(size = 12), panel.grid.major = element_blank(),
          panel.grid.minor = element_blank(), plot.title = element_text(size = 12),
          plot.margin = unit(c(top, right, 0.5, left), "cm")
    ) +
    geom_boxplot() +
    coord_flip() +
    scale_fill_manual("", values = sampleCols) +
    ggtitle(paste0("Pvalue=", bPValue))

  legend <- g_legend(p + theme(legend.position = 'left', legend.text = element_text(size = 8),
                               legend.margin = margin(t = -2, unit = 'cm'))
  )

  gA <- ggplotGrob(p)
  gB <- ggplotGrob(bBoxPlot)
  maxWidth = grid::unit.pmax(gA$widths[2:5], gB$widths[2:5])
  gA$widths[2:5] <- as.list(maxWidth)
  gB$widths[2:5] <- as.list(maxWidth)

  grid.arrange(gA, rBoxPlot, gB, legend, ncol = 2, widths = c(0.75, 0.25), heights = c(0.75, 0.25), newpage = F)

  dev.off()
}


